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ABSTRACT 


This thesis describes the analysis, method and computer 
programs used to obtain the Fast Fourier Transformation (FFT) 
of an electroencephalogram (EEG) using a small laboratory 
computer like the, PDP Lab 8/E. The EEG power spectrum was 
then computed from this transformation. The information 
contained in this thesis is intended to enable the user to 
compute the Fourier coefficients of a set of data points or 


compute the power spectrum of a real waveform such as the 


EEG. 
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I. INTRODUCTION 


As. cerebral cortex cells discharge electrically, a small 
potential is produced; when this potential is averaged with 
the ee discharging cells it can be sampled with electrodes 
applied to the scalp. This weak voltage is amplified 100,000 
times by a low noise, high gain amplifier. For simplicity 
and economy, a paper and pen writing oscillograph is used to 
record the changing potential.. This record of cortical po- 
tentials is called an electroencephalograph (EEG). 

These waves are very complicated and analyzing them has 
been, a. problem since, their.discovery.. Due to the wave-like 
nature, EEG frequency or spectrum analysis has been the 
subject of much research over the years. Large, high speed 
computers are being used to aid this work today. The pur- 
pose of this study was to develop fast Fourier techniques 
for spectrum analysis of EEGs using a small laboratory com- 
puter like the PDP Lab 8/E computer system. With the com- 
puter located in the laboratory, near subject and equipment, 
it was hoped that some kind of real time frequency analysis 


could be accomplished. 
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II. BACKGROUND 


A. PHYSIOLOGY 

The small potential measured in the EEG originates in 
the outer .1 inch of the brain called the grey matter. The 
working cells of the grey matter are the neurons, numbering 
perhaps to ten billion in all. Very little is known about 
how these cells operate to provide the many functions that 
the brain accomplishes. One well known and supported theory 
by John C. Eccles explains how the EEG waveforms are pro- 
duced [Ref. 10]. Each neuron is assumed to be an independent 
functioning unit, receiving stimuli from outside sensory 
devices and other neurons through dendrite fibers. When a 
particular number of stimuli are received, the neuron fires 
a pulse of electricity to other neurons via axon fibers. 
Beecerving reverse polarity signals will cancel the effect of 
incoming stimuli and prevent the neuron from firing. In 
this way each neuron may be connected to as many as 100 
other neurons. 

Usually many thousands of neurons are stimulated simul- 
taneously. As they reach their excitation points and fire, 
other thousands of neurons receive these new stimuli. Thus 
the collection of stimuli moves through the brain step by 
step. Each stage in the path may be only a thousandth of 
an inch or it may be three or four inches in another part 
of the brain. After firing, each neuron has a deactivated 


period of several seconds, in which it cannot be stimulated 


7 













ar eyo en Fgia oan ens it bekensem m 


fio yorg siz baifes — st 


= 


om 


gr i . If it OF 
x7 
je ; Ls 2 ¥ 
sid = 7; obivorTd 
Pte j ¥ srt 


, an ee 
foimen doad 

griwisset ,1 

Say. sao wen 

va Elumize to 13d 

4a 

5 oF Yatolrtsoss 

; nee 

mie ~rivatog server 

rong bre it soe 8 


yipewouds yee 


a 


oJ yr ae ah 


ytuen io coer 


= hy og 


fo respond. This delay, or refractory period as it is called, 
produces a wave-like quality to the stimuli collection as it 
moves along. Eccles believes that this explains the wave- 
like nature of cortical potentials as they filter out through 
the skull and scalp and are recorded on the EEG. 

Another important aspect of these neural pathways is that 
they frequently double back on each other. This causes a 
circular path to form, which provides a self-exciting fea- 
ture to the brain. These closed loops, with fixed distance 
and refractory periods, may account for the predominance of 
certain frequency waves. The prospect of neurons having 
feedback loops opens up a whole new level of complexity to 
the brain. In this way each neuron could be viewed as a 
very flexible amplifier, each one operating differently, 


depending on the kind of feed-back it is receiving. 


B. HISTORY 

A short history of the cortical potential must begin 
with Richard Canton, a British physiologist, who reported 
in 1875, ''Feeble currents of varying direction pass through 
the multiplier when electrodes are placed on two points of 
the external surface (of the brain), or one electrode on the 
grey matter, and one on the surface of the skull." At the 
time, Canton was examining cat, monkey and rabbit brains 
with a sensitive galvanometer using optical magnification. 
It is remarkable that he could distinguish such small po- 


tentials without the use of electronic amplification, which 
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wasn't developed for another fifty years. In 1876 the 
Russian, Danilevskey, discovered a change in the average 
potential of the cortex in response to acoustic stimuli. 
Only recently has this interesting aspect of brain activity 
been reinvestigated. In 1914 Cybulski recorded an apileptic 
seizure in a dog. The early attempts at recording cortical 
potential used photography, however, this was expensive and 
tedious. Pen recorders were not available for this use 
until 1940 [Ref. 8]. 

Recording cortical potentials in man lagged behind their 
demonstration in animals for many years. In man, both elec- 
trodes were attached to the scalp, causing the potentials 
to be greatly attenuated through the skull and scalp. The 
first recordings from the human brain were made in 1924 by 
Hans Berger in Germany. He worked and published much, ex- 
ploring many aspects of the use of the electroencephalogram. 
Some of these included physiological, neurological, psychia- 
tric and psychological applications. He confirmed Canton's 
findings in animals and established that these waveforms 
change with age, with sensory stimuli and with chemical 
changes in the body. He was the first to record a major 
epileptic seizure in man. Due to the earlier work by 
Kaufman and Cybulski, it became known that epilepsy in 
animals caused abnormal waveforms. 

After Berger's first publication in 1929, many labora- 
tories took up this. work and knowledge in the field grew 


rapidly. After substantiating Berger's claims, the various 
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brain centers began to specialize in different research 
areas: Gibbs and Lennox in epilepsy, the Davises in normal 
subjects, Adrian, Hoagland and Jasper in the physiological 
mechanisms and Walter in brain lesions. The growth of 
electronencephalography has really paralleled the develop- 
ment of the machinery for recording the waves. In the 1930's 
the electronic amplifier gave a tremendous boost to the 
field. In 1948 the transistor allowed the development of 
highly portable equipment. Modern instruments allow 16 
channels to be monitored simultaneously, each of which has 
linear response and accurately reproduces frequencies from 
1 to 500 Hz. A more complete history of electroencephalo- 
graphy may be found in Braizier [Ref. 4]. 

One of the many things that Hans Berger discovered about 
the human electroencephalograph was that the most normal 
pattern in adults, with eyes closed, was a frequency between 
8 and 13 Hz which he called alpha rhythm. These and other 
rhythms seem to originate in the cortex, however, they cease 
if the cortex is surgically disconnected from the thalamus. 
Alpha rhythm is most strong when the subject is relaxed 
With eyes closed. Opening the eyes in normal light signifi- 
cantly reduces alpha rhythm and similarly, opening the eyes 
in total darkness eagaes altéemporary reduction. There is 
some evidence [Ref. 18] to support a positive correlation 
between high amplitude alpha rhythms and a passive or in- 
troverted personality structure as revealed by psychoanal- 
ysis. Similarly, low amplitude alpha rhythms were found 
prevalent in extroverts [Ref. 18]. 


10 





i ee 


Das ¢ 


i: 
ie) ‘ 
is i, 


— 
i” 
sie 


 anevee ok aan ah ae rons a a9 


2/006 af? nf <eevar oat: peiilirosen tot ytsatdoamt 





























i aeeonacnen ne yequad baw bncigeod - & c 
a o diworg oT eto hdr ab ths ot Bam 3 
-qotevab odt befottetig “btnvy eed yigergo Lang 





903 ot teoed awobasti1? « eveg trolbilqas ofa 

to Gnsemqolerob ott bowolia r1dde2 Penats ait? anen t 
at wolls etosmusteni nrobol! . tnomgtupe stds 

‘eed AOdR lo fone .yfawoomitingie hoerot tom ad's 
work estonaupstt essvhotier eletatu 23% tas 
-oladigannoeorsosio. Yo yroselh ssetQme> .oTor A att 
e; DA ,¥o8). vabstoerd wh Bauet os , 
unde hte ovooe tS Yopwed: ete sone agente? a ee ae 
tunton txom on} tniy sow iqatgeledasonsni sole em at 
reevied Yoreupst? n ew \beecls coy taiw ,20iebe mie 
yeito haa weet .colavde edate Sbalian 44 ied iG “i: 


; J 
xXoltgs otfent signigity® 62 age 


S265 20 rOVSWOl 
sweated) etfs wort be tosaneu4®> Vhigokyuwe Se rood 
Wexefor 2k tosttie sao noilw grovee Ta0n 22 WE 7 
note wigqit 4 i pier -beeol> 28Ne 
xo Sa? gatnego Cig nix sigtx esduber re 
ston i touh: : une eapndish fates 2 


2) . deh) oeerebive ae 
} Yi Reon P.eng eri iy Ps Lents hota nosweed - 


uy 


Delta activity is between 1 and 3.5 Hz and is generally 
associated with sleep. In 1954 Bremer achieved a strong 
delta state by disconnecting all the sensory pathways. This 
is logical since sleep is most easily induced when sensory 
stimulation is reduced to a minimum. It must be emphasized 
however, that even though no light, sound, etc. is sensed, 
one does ‘not 'sieep unlessthe is "tired," i.e., delta state 
and sleep are not the same. Theta activity, in the 4 to 7 
Hz range, is usually associated with normal children [Ref. 4]. 
The percentage of theta decreases significantly with age in 
normal subjects. Beta activity, with a frequency of 14 to 
30 Hz, is generally associated with increased anxiety levels. 
The incidence of beta is increased in subjects under-going 
alcohol withdrawai. Beta activity can easily be blocked in 
the motor center associated with a particular movement [Ref. 
4]. 

The four basic waves, alpha, delta, theta and beta, span 
the frequency spectrum of cortical potentials. Another pat- 
tern of interest is the mu rhythm, which is characterized 
by a sharp peak in negative half-cycle, rounded positive 
half-cycle and 7 to 11 Hz. Mu rhythm is reduced not only 
by overt movement of a limb but also in some cases when 
the subject is contemplating movement [Ref. 4, 8]. 

Except for beta waves, all the waveforms previously 
described tend to be reduced by some kind of stimulation. 

In contrast, lambda waves, a saw-toothed waveform, are best 


seen when the subject is actively looking at something of 
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interest to him. This waveform usually occurs when the 
subject is scanning a complex pattern or picture. Paradox- 
ically, random waves of the same polarity and shape occur 
in some subjects during light sleep. A thorough analysis 
of the various waveforms in psychologically normal and 
abnormal subjects may be seen in Thorp and Heyman [Ref. 24] 
and Kiloh and Osselton [Ref. 21]. 

Interpretation of the EEG by visual inspection of the 
record depends in a large part of recognizing common pat- 
terns and their association with psychologically normal or 
pathological conditions. Unfortunately these visual methods 
are unable to differentiate between clinical conditions 
having similar EEG's. Several different methods for anal- 
Melee the Eau dre being tested. "Une ot the’ first attempts 
was to get a quantitative grasp on the amplitude of the 
waves. Walter [Ref. 4] joined the peaks and troughs with 
two lines and then measured the average distance between 
the cee With the advent of amplifiers, electronic inte- 
gration solved this problem. 

In looking for better ways to specify cortical poten- 
tials, many forms of frequency analysis have been tried. 

As early -as 1932, Dietsch suggested that these potentials 
might be composed of sine waves and he outlined the procedure 
whereby Fourier analysis could be used to determine the fre- 
quency components present [Ref. 9]. 

The first frequency analysis of an EEG was accomplished 
at the Cincinnati College of Medicine in 1944, where the fre- 
quencies were determined by manually counting the number of 
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complete waves in each second [Ref. 12]. In the late for- 
ties and fifties, much work was stimulated by Gibbs and 
Grass [Ref. 13] when wide-band electronic filters were used 
te separate the different frequencies [Ref. 19, 22]. Filter 
devices are commonly used in this way, even today [Ref. 25]. 
An excellent review of all work on EEG analysis techniques 
up through the fifties may be found in Reference 5. 

Spectrum analysis of physiological data using the dis- 
crete Fourier transform was accomplished by Halberg in 1961 
{[Ref. 14, 23]. Walter was first to use the discrete trans- 
f£OFMs ony, cCortHeal spotentials;in 1963; [Ref.s, 26, 28]. ,Most of 
the quantitative frequency work being done on EEG's today 
is with ie aid of large computers in which the discrete 
Fourier transform, provides, the.conversion, {Ref.,11,.15, 17]. 
ehetr v.15: Navy is also entering this o research) area, rat. the 
Navy Medical Neuropsychiatric Research Unit in San Diego, 
California [Ref. 15]. In 1973 the Fast Fourier Transforn, 

a faster version of the discrete Fourier transform, came 


into use on the EEG [Ref. 2]. 


ES 


















ii a ie ae tai nee 
; mb lage oT 1 Me 


pbs - ail bi Py y ih ph 


rats 


a I Bi sig | 
i iit. monk 
" nh aon oe = »: 
seitoloe’, sieges one ne. tno ‘tat ey 
7 sotarol of at bawg? od Yom aoe 
ath et guteu stab teobgoloteda to sia iam 
toes tf gisdi an xd Sorludiqneass esw voted 
“aks? aperseth only say or texkt ame eRe 
1a 220M . {8S OS: . Per) FORE Bh, thelsneteg, 198 
ynbos acre 19 @fob ghied drow Yoaeupet? syits 
ore sath ada dobilw- sg eros uqmNeD agtef to bia 
ETL . ~tt-- EE. 2. BoLesoveos at: -8ob ive g-mEOT ERO El 
od? te £918 tovesesr 22nd ankredas oats ki aha 
ogni gee ai Jisl detest 9iveehtoyeqaneen a bbs 
ivotentba? :ofayel t2g% sdf OVOP gE. ioe ee 1 
subo ,mwrobernst) tainted siegoete Sas ae noiaroy, 19m 

[2 .2oR]. S82 of2 ae 


III. FOURIER ANALYSIS 


A. FOURIER TRANSFORM 
1. Continuous 
An introduction to Fourier analysis must begin with 


the general form of the Fourier series. It is 


E(t) = Aj + ,2, A, cosnwt + .2, B, sinnot (1) 


where w = 21/T and T is the period. Any waveform can be 


represented by determining the coefficients Ay? A_ and Bo 


n 
and summing all the terms of this series. The coefficients 
give the relative contribution of each frequency that is 
present in the signal. These are hermonics or integer 
multiples of the fundamental frequency 1/T. Fourier de- 


veloped a set of equations from (1) and solved them simul- 


taneously. The result 


A= ‘i £(t)dt (2) 

A. = 2 A f(t) cosnwtdt (3) 
n Il, f 

ee £(t) sinnwtdt (4) 
n alt sf 
oO 


provides a method for finding the harmonics when the function 
is known. 
2. Discrete 
Frequently, however, the complete function is not 


known and only a discrete number N of samples over a finite 
14 
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period T is known. In this situation, an approximation to 


the continuous function may be written 


n n . 
f(t) = Ag + 2 el Ae cosmwt + 2nz1 BF sinmwt (5) 


where n = N/2 for N even and n = N+1/2 for N odd. The ex- 


pressions for the coefficients now become 





ee Pee! 27mt 

An = po Fit ) Cos ( Ai ae ) (6) 
Aes : 27mt . 

Bh ae Le. HC He sinovwG NV iM (7) 


These functions are easily computed from the set of 
sample data points; however, the time required is sometimes 
prohibitive. It will be noted that this procedure requires 
N? trigonometric operations, since each of N coefficients is 
computed from N terms of a series. These computations, in 
addition to numerous other arithmetic operations, require 
considerable computer time, especially on small laboratory 


computers like the PDP Lab 8/E. 


B. PROBLEM AREAS 

There are three problems commonly encountered when using 
the discrete Fourier transform. These are: aliasing, leak- 
age and the picket-fence effect. Aliasing occurs when a 
higher than expected frequency is sampled, but appears in 
the spectrum as a lower frequency. Suppose an unknown sig- 
nal is being sampled five times at a rate of 100 Hz. If 


an 80 Hz component is present, it will appear in the spectrum 
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as 20 Hz and there is no way of distinguishing this pseudo- 
frequency from the real one. Nyquist proposed the solution, 
in that, the sampling rate should be twice the highest fre- 
quency present [Ref. 16]. This can'be accomplished in an 
unknown signal by passing it through a lowpass filter set 

at half the sampling rate. A lowpass filter was used to 
prevent aliasing in this study. 

The problem of leakage occurs in all Fourier analysis 
using finite sampling records. The record, or the total 
length of the sample, can be viewed as a window. When this 
window is transformed, without the sample data present, the 
resulting spectrum contains not only the fundamental fre- 
quency but also spurious peaks called sidelobes. When the 
window with its sampie is transformed, these sidelobes tend 
to leak into the spectrum. Because of the relatively small 
sidelobes, in the frequency range of interest, leakage was 
not a problem in this study. 

The picket-fence effect is a problem when a pure sinu- 
Soid occurs between two adjacent harmonics. If it is ex- 
actly halfway between computed harmonics, the signal is 
reduced to .637 of the original in both harmonics. The 
worst situation occurs when there is a frequency between 
each pair of adjacent harmonics. In this case, an erroneous 
ripple, or picket-fence, appears in the spectrum. Pure 
Sinusoids were assumed not to be present in EEG's and thus 
the picket-fence effect was ignored. Bergland discusses 
these problems as well as the mathematics associated with 
them in Referénce 3. 
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IV. FAST FOURIER ANALYSIS 


A. FAST FOURIER TRANSFORM 

In 1965 James Cooley and John Tukey published an impor- 
tant paper that described a much faster method to compute 
the discrete Fourier transformation, called the Fast Fourier 
Transform (FFT) [Ref. 7]. The procedure adopted for use on 
high speed computers has become known as the Cooley-Tukey 
Algorithm. A general outline of the proof of this algorithm 
follows. 

A time series Xy containing N even samples, can be par- 
titioned into two functions : and Ze each containing N/2 
poinces., "Let 2a be composed of the even-numbered points 
(X, 9X29+-+ sXyy) and Z.. the odd-numbered points (XqoXze--Xy_yy2)- 
Equations (6) and (7) may be combined using complex number 


notation, so that 


Ke wit (8) 


where W = eon aga 49 = ¥ =i... Since = and Z.. are a Set of 
N/2 points, each has a discrete Fourier transform whose 


coefficients are 


N/2-1 ae 
Bh h rZ0 Y, a (9) 
N/2-1 
i m(2r+1) 
Ch = rz0 ae W . (10) 
N 


with m = 0,1 1. Combining equations (9) and (10), to 


Tet a7 
get odd and even points, gives 
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The details of this proof will be found in Reference 6. 

If Xx. is an even power of two, then ss and oe also have 
an even power of two points and each can be partitioned again 
and again until only one point remains in each set. The 
Fourier transform of a single point is just that point again; 
the entire transformation can be accomplished using equations 
(11) and (12) repeatedly at each partition. The requirement 
that the number of points N be a power of two assures that 
after log, N partitions there will be only one point in each 
set. In this way, the entire set of points can be translated 
without carrying out the manipulations of equations (6) and 
(7). Since the FFT is derived from the discrete Fourier 


transform, the FFT maintains the same precision as before. 


B. FLOWPATH 

An example of the FFT for N = 8 will be illustrative. 
Figure 1 shows the flow diagram and the configuration after 
Bae tiist partition. | Using Tukey's notation, each dot rep- 
resents a variable and the factors entering the dot are 


additive. When a multiplication is required before adding, 
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Fast Fourier Transform, Second Partition, Flow Graph. 
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the quantity will be beside the factor to be multiplied. 

The set of ie for j = 0,...,7 is computed from equations 

oa) ana: (12). The quantity near the bottom right dot is 

AW = B, + W’C.. Figure 2 shows all quantities after two 

partitions and Figure 3 shows the complete transformation 

broken down into only multiplication and addition operations. 
It may be noted in Figure 3 that the data points must be 

Fearrangped according. to the bit reversal rule before the 

algorithm can be carried out. Bit reversal is accomplished 

by representing the data point position number in binary 

and reversing: the bits about the center. Position (3) is 

mwreceen (0 Wi) yand is reversed to (1 7 0), which. is. (6). 

Thus (0) stays a (0), (1) becomes a (4) and so on. More 


information about FFT may be found in-References 3 and 6. 


C. PROGRAM 
The computer program FAST was written to accomplish the 
FFT on the PDP Lab 8/E computer system. A complete listing 
of FAST may be found in Appendix A. The program is written 
in BASIC, one of the commonly used languages on the PDP Lab 
8/E. FAST receives data points as inputs and outputs the 
Fourier coefficients. 
In order to operate the FFT, read FAST into core using 
normal methods, and een ae RUN. The program will type 
ENTER DATA POINTS 
? 
and go into a wait loop. Type in 64 data points in their 


proper. order. The number of data points may be changed to 
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any power of two, but cannot exceed 64 due to core limita- 
tions, which is 8,000 words in the PDP Lab 8/E. FAST will 
translate the data points and type the Fourier coefficients. 
The time required to convert 64 points is approximately 25 
seconds. The power associated with each harmonic can be 
obtained by squaring and adding the coefficients at each 


frequency. 


D. EXAMPLE 

As an example to illustrate the use of FAST, the square 
wave shown in Figure 4 was translated and the results com- 
pared with the coefficients from the continuous Fourier 
transform... Tabie I shows the result of a FAST translation 
of the square wave along with the coefficients. To deter- 
mine the A coefficients using the continuous transform, 
equations (2) and (3) are integrated from -T/2 to T/2. 
Since £(t) is an odd function, the AL terms; arewall-zere. 
Using equation (4), however, and integrating over the same 
range, it is easily shown that 


Bp = 22. 2Cosnt (13) 
n 7 Sait 


To compare with FAST, the same 32 harmonics were computed 
using equation (13)(see Table II). The answers from FAST 
are within two decimal places of the answers computed from 


the continuous Fourier transform. 
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Figure 4. Square Wave. 
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Figure 5. EEG Wave and Timing Pulse. 
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att bas ovat OSa ae origil 


HARMONIC A(N) B(N) 


0 0 0 

1 0 P27 2217 

2 0 0 

2 0 -4213407 
4 07 0 

5 (0 - 2495139 
6 0 0 

7 0 - 1746758 
8 0 0 

= 0 sfSZ1451 
10 0 0 

11 0 -1042749 
1.2 0 0 

a5 0 -08427145 
14 0 0 

1S nee - 06895806 
16 0 0 

iW 0 .05664668 
18 0 0 

19 0 - 04635313 
20 0 0 

as 0 -03746103 
22 0 0 

23 0 -02956024 
24 0 0 

ye 0 02236283 
26 0 0 

oa 0 sG1565 959 
28 0 0 

29 0 Oi.2 709S8E=035 
30 0 0 

31 0 S.070205E-03 
READY. 


Table I. FAST Coefficients Translated from the Square Wave 
of Figure 4. 


24 





ecreebs. 
szcaett, 
reases. 
eprsane, 
ansrc88e, 
avseenan, 


eGopedad. 
f 


secsocaocecorescooeooessoee seesoseseey $ 


~ 


met i 
i” ae : 
TOM ie Mien ri 
; nee). 18y 
tar. oh 
a) ee 

ve a 
me 





'a¥e0) TEA? 
bh osug it te 


RUN 
HARMONIC 


WOONKAMNRPWNH © 


Table II. 


B(N) 


-507059E+37 
~2tacel 

LI eLole=09 
~4244135 

.o i590 50-09 
- 2546481 
-581064E-09 
LG LEILS 
t69/798E-09 
~1414712 
-486382E-10 
tLe 7eoe 
-905318E-10 
-09794159 
ati ooo.e= LO 
.08488271 
«928989E-10 
-07489651 

ney a Uae WA Ba) 
VO70LZO7 
—L45191E- 10 
.06063051 
~JLL9928-10 
~ 05955629 
aoa20 598 = 10 
~ 05092963 
.648608E-10 
-04715706 
soa) 99 5E— 10 
-04390485 

a hi2 L27E= 10 
.04107228 


Continuous Fourier Transform of BL 
for Square Wave of Figure 4. 
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V. SPECTRUM ANALYSIS 


A. SPECTRUM 

In order to do frequency analysis of EEG waveforms, a 
computer program called SPECTRUM was written which incor- 
porated FAST (Appendix A). The problems encountered were 
many. First, a relatively noise-free signal had to be pre- 
sented to the PDP Lab 8/E analog-to-digital (A-D) converter. 
Then the signal had to be sampled at the right rate and 
duration. Finally, the data had to be delivered to FAST in 
a compatible format for translation. The power spectrum was 
chosen as the output format since it has more intuitive ap- 
peal than some other measure such as voltage or current. 
According to Parseval's theorem, the power at each frequency 
can be computed by squaring the two coefficients and adding 
them together. This was done and the power spectrum is dis- 
played on the teletype as percent of total power at each 


harmonic. 


B. APPARATUS 

The Beckman type RM Dynograph Recorder was used to pick- 
up and amplify the cortical potentials. Using sodium chlo- 
fade paste, silver electrodes were applied to the subject's 
scalp at positions Cy and T6: A standard Beckman A-C 
coupler was used with the TIME CONSTANT switch set to .03 
and HIGH FREQ switch on 3. Depending on the signal strength, 


the preamplifier was set to either 1 mv/cm or .5 mv/cm and 
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the power amplifier set for a multiplication of .02. This 
adjustment was made in order to maintain a relatively con- 
seat Signal level from the power amplifier. The FILTER 
switch was set on 3, which is a 60 Hz lowpass filter. 

The EEG was recorded on one track of a Hewlett-Packard, 
four track, instrumentation recorder on 1.0 mil magnetic 
tape at 3-3/4 in/s tape speed. A seven ms timing pulse was 
recorded on another track from an Ortec pulse generator. 
The timing pulse made it possible for the computer to find 
a particular place on the tape. The two tracks were then 
connected to the A-D converter of the PDP Lab 8/E computer. 
As an alternative hookup, the A-D converter can be conneeea 
directly to the Beckman and sample the EEG as the subject 


produces it. 


C. SAMPLING RATE 

With the proper signal,on. channel. 1. of. the A-D converter, 
SPECTRUM samples and translates it to the power spectrum. 
According to Nyquist's rule, the sampling rate must be twice 
the maximum frequency. The maximum frequency passed by the 
Beckman filter is 60 Hz, so the sampling rate must be at 
beast.) 20 Hz... lm order, to. give. Jittlesextna margin, 128, Hz 
was chosen so that the harmonics would come out in even 
multiples of two. FAST can only convert 64 points and sam- 
pling at 128 Hz, the period (T) of the sample is % second. 
The reciprocal of the period (1/T) determines the increment 
of the harmonic. In this case 1/T = 2, thus each harmonic 


is exactly 2 Hz apart. 
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D. EXAMPLE 

After loading SPECTRUM and typing RUN, the program sam- 
ples, converts the data, and types the power spectrum. 
Table III shows an example of one SPECTRUM run, taken from 
a subject who was relaxing with his eyes closed. The pro- 
gram starts sampling when it receives a numbered timing 
pulse on channel three. The EEG waveform from which this 
example was taken, along with its timing pulse, is shown 
werigure 5. SPECTRUM. uses the timing pulses to locate the 


desired sample on the tape. 


E. EXPERIMENT 

By choosing the proper segment with the timing marks and 
replaying the tape 20 times, SPECTRUM was used to scan a 
.length of tape arbitrarily set for 10 seconds and produce 
the power spectrum for each half second record. The results 
of this experiment were reduced to the graph of Table IV. 
It plots frequency against time, with each number represent- 
ing the percentage of total power for a given frequency and 
time. If the table is viewed as a three dimensional graph, 
where percentage power is shown coming out of the page, it 
ZS Gasy to visualize a ridge running up the 10 Hz line. The 
subject ae again resting with eyes closed, so the prepon- 
derance of alpha frequency (8, 10, 12 Hz) was understandable. 
The frequency also appeared to come in bursts. It was first 
on for about three seconds, then reduced for one second, then 


back on for a second, and so on for the length of the tape. 
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RUN 


-FREQUENCY ; POWER 


HZ % 

0 16.2738 
1 1.157298 
2 .07388068 
3 4.684211 
4 5.177504 
5 3.729045 
6 49.26551 
7 5.911107 
8 2.511872 
9 3.593473 
10 - 7096686 
11 -6507585 
12 -9833005 
13 .9566361 
14 -6030759 
15 -538035 
16 -2399192 
ily -2305483 
18, « - 3325844 
19 .2973473 
20 .2739161 
7A .1222699 
22 . 2039538 
235 . 178039 
24 . 1704812 
25 .1763688 
26 - 1824439 
Di . 1651099 
28 - 1600951 
29 .1526896 
30 . 1463384 
31 .148733 
READY 


Table III. A Typical SPECTRUM Run from a Half Second Sample 
of an EEG; 
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17 

9 
3 

8 
8 
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13 
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FREQUENCY (Hz) 


Power Spectrum derived from an EEG Over a Period of 
Ten Seconds. “Each Number Represents the Percentage 
Power at each Frequency for the Half Second Sampled. 
Zeros are’ omitted, 
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Sometimes this alpha burst effect can be detected right on 
the EEG.” ’ 

There were two places where strong delta frequency was 
present. Prepresented on the graph by ‘the 2 Hz line, ‘delta 
was present in the beginning and again around the five 
second point. This normally indicates a very light stage 
of sleep, where the subject is drifting in and out of ‘con- 
sciousness. Under questioning the subject indicated that 
he could have been asleep several times during the recording 
session. The graph also shows the rapidity with which fre- 
Guency changes [Ref. Z; 11, 15 and °27]'. ° One ‘example of this 
is from the first half second to the one second point, where 
the 10 Hz power changed from 3% to 39%. Rather than changing 
smoothly, it was as if one frequency stopped and another 


started. 


F. NOISE 

Either because of a noisy amplifier or possibly a dirty 
tape, it may become necessary at some time to analyze a 
noisy signal. This will become known when repeated spectra 
from the same point show a variance. In order to reduce 
this variance, a method of averaging was devised. Instead 

of Sanpling at 128 Hz, a new rate of 1280 Hz was. used and 

the sample means were computed every 10 points. The sample 
period was kept at & second, so a total of 64 sample means 
were computed each run. This reduced the variance by a fac- 
tor of Y10 and providing the signal’ was not too noisy, the 


results were satisfactory. Thus to obtain more information 
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from the noisy EEG, it may become necessary to sample more 
frequently than according to Nyquist's law. 

Increasing the sampling rate is not without its short- 
comings. Firstly, it takes longer to process the added data. 
Changing to a 10 point average increases the computing time 
Eon 25 to 40 seconds. Secondly, part of the data must be 
stored while the averages are being computed. The new pro- 
gram requires about 200 more storage positions, which is at 
the maximum capability of the PDP Lab 8/E with 8,000 word 
memory. Most spectral analysis work being done on EEG's 


today use a 200 Hz sampling rate [Ref. 2, 15]. 
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VI. CONCLUSIONS AND RECOMMENDATIONS 


A real time analyzer for the EEG was one of the goals 
of this study. This has only partly been realized since 
SPECTRUM requires 25 seconds to convert a half second sam- 
pre. Since the EEG frequency pattern changes So. rapidly, it 
is necessary to convert more than a half second; perhaps as 
mech as 10 seconds in order to get a good representation of 
the changing nature of the wave. This means. that a much 
larger computer is needed to do real time spectral analysis 
bE EEG's. 

As was shown using off-line techniques, SPECTRUM is a 
useful tool for spectral work. Using an instrumentation 
recorder in conjunction with some kind of pulsing device 
for timing, many different forms of frequency work may be 
accomplished such as autocorrelation, cross spectral analysis 
and digital filtering. SPECTRUM can also be used to compute 
power spectra from other kinds of continuous waveforms such 
as electrocardiograms and electromyograms. 

The author believes that this area of study is an impor- 
tant one and should be extended in the future; however, 
there are limitations. The major limitation to this system 
resides in the lack of adequate data storage, both in the 
small core size, which is about 4,000 words with BASIC in 
core, and in the lack of a random access device such as a 
disk pack. Of these two, the disk pack would be most useful 


for spectrum analysis work since large amounts of data from 
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the A-D converter could be stored rapidly and processed as 
required. The flexibility of the entire system would be 

considerably enhanced because the disk pack could also be 
used to break up large programs which are now too big for 


the computer and then process them in a segmented fashion. 
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APPENDIX A 


SPECTRUM 
E N=64 
2 N1=63 


5 DIM F(1,63),W(1,31),T(1) 
8 SET RATE 6,558 

10 DIM P1(300) 

15 REAL TIME P1,1,0,900 
18 IF ADC(3)<=.9GO TO 18 
20 ACCEPT 

30 FOR I=0 TO 63 

35 FOR J=1 TO 14 

40 F(0,1)=F(0,1I)+ADB(0) 
45 NEXT J 

50 NEXT I 

58 PRINT 

100 FOR I=0 TO 31 

110 READ F(1,1) 

115 F(1,63-1)=63-F(1,I) 
120 NEXT I 

158 FOR I=0 TO N1 

162 IF F(i,i)=I THEN 169 
163 X=F(0,I) 

164 J=F(1,I) 

165 F(0,1)=F(0,J) 

166 F(0,J)=x 

167 F(1,1)=F(1,J) 

168 F(1,J)=J 

169 NEXT- I 

194 FOR I=0 TO N1 

195 F(1,1)=0 

196 NEXT I 

210 W(0,0)=1 

220 W(0,1)=COS(6.283185/N) 
230 W(1,1)=-SIN(6.283185/N) 
240 IF N<=4 THEN 295 

250 FOR I=2 TO N/2-1 

260 I1=I-1 

270 W(0,1)=W(0,1)*W(0,11)-W(1,1)*W(1,11) 
- 280 W(1,I1)=W(1,1)*W(0,11)+W(0,1)*W(1,I1) 
290 NEXT I 

310 J=1 

320 FOR M=1 TO N 

330 J=J*2 

340 IF J>=N THEN 370 

350 NEXT M 

370 FOR K=1 TO M 

380 J2=2+(K-1) 

390 I2=N/(J2*2) 
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400 FOR I=0 TO 12-1 
401 T(0)=F(0,Y) 

410 FOR J=0 TO J2-1 

500 X=N*I/I2+J 

510 Y=X+J2 

520 Z=J*I2 

530 T(0)=W(0,Z)*F(0,Y)-W(1,Z)*F(1,Y) 
540 T(1)=W(1,Z)*F(0,Y)+W(0,Z)*F(1,Y) 
550 F(0,Y)=F(0,X)-T(0) 

560 F(1,Y)=F(1,X)-T(1) 

570 F(0,X)=F(0,X)+T(0) 

580 F(1,X)=F(1,X)+T(1) 

590 NEXT J 

600 NEXT I 

610 NEXT K 

690 A=0 

700 FOR I=0 TO N/2-1 

710 F(0,1)=(F(0,1)+2)+(F(1,1)+2) 

712 A=A+F(0,1) 

720 NEXT I 

728 PRINT "FREQUENCY POWER" 
730 FOR I=0 TO N/2-1 

732 A2=F(0,1)/A 

750 PRINT I*2,A2 

760 NEXT I 

800 FOR I=0 TO Nl 

810 F(0,1)=0 

820 NEXT I 

900 DATA 0,32,16,48,8,40,24,56,4,36,20,52,12,44,28,60,2 
O0i DATA 34, 18,50,10,42,26,58,6,38,22,54,14,46, 30,62,1 
998 RESTORE 

999 GO TO 1 


1 N=64 

2 N1=63 

5 DIM F(1,63),W(1,31),T(1) 
30 FOR I=0 TO N1 

35 INPUT F(0,I) 


49 F(1,1)= 
50 NEXT I 
58 PRINT 


100 FOR I=0 TO 63 

. 110 READ F(1,I) 

120 NEXT I 

130 FOR I=0 TO N1 

140 11=1*64/N 

150 F(1,1I)=F(1,I1) 

155 NEXT I 

158 FOR I=0 TO N1 

162 IF F(1,I)=I THEN 169 
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163 
164 
165 
166 
167 
168 
169 
194 
195 
196 
210 
220 
230 
240 
250 
260 
270 
280 
290 
310 
320 
330 
340 
350 
370 
380 
390 
400 
401 
410 
500 
510 
520 
530 
540 
550 
560 
570 
580 
590 
600 
610 
698 
699 
700 
710 
720 
730 
736 
738 
739 
740 
760 
761 


X=F(0,I) 

J=F(1,1) 

F(0,1)=F(0,J) 

F(0,J)=X 

F(1,I)=F(1,J) 

F(1,J)=J 

NEXT I 

FOR I=0 TO N1 

F(1,1I)=0 

NEXT I 

W(0,0)=1 
W(0,1)=COS(6.283185/N) 
W(1,1)=-SIN(6.283185/N) 

IF N<=4 THEN 295 

FOR I=2 TO N/2-1 

I1=I-1 
W(0,1)=W(0,1)*W(0,I1)-W(1,1)*W(1,I1) 
W(1,1)=W(1,1)*W(0,11)+W(0,1)*W(1, 11) 
NEXT I 

J=1 

FOR M=1 TO N 

J=J*2 

IF J>=N THEN 370 

NEXT M 

FOR K=1 TO M 

J2=2+ (K-1) 

I2=N/ (J52*2) 

FOR I=0 TO I2-1 

T(0)=F(0,Y) 

FOR J=0 TO J2-1 

X=N*I/I2+J 

Y=X+J2 

Z=J*12 
T(0)=W(0,Z)*F(0,Y)-W(1,Z)*F(1,Y) 
T(1)=W(1,Z)*F(0,Y)+W(0,Z)*F(1,Y) 
F(0,Y)=F(0,X)-T(0) 
E(1,Y)=FC1 ,X)-T(1) 
F(0,X)=F(0,X)+T(0) 
F(1,X)=F(1,X)+T(1) 

NEXT J 

NEXT I 

NEXT K 

A=0 

PRINT "HARMONIC A(N) B(N)" 
FOR I=0 TO N/2-1 

PRINT I,F(0,1)/32,F(1,1)/32 
NEXT I 

FOR I=0 TO N/2-1 

V1=2*I/N 

V2=F(0,1)/A 

IF V2<1.000000E-03GO TO 760 
PRINT I*2;V2 

NEXT I 

STOP 
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800 FOR I=0 TO N1 
810 F(0,1I)=0 

820 NEXT I 

B00 DATA 0,32,16,48,8,40,24,56,4,36,20,52,12,44.28,60,2 
BOL DATA .34,18,50 ,10,42,.76.58.6,38,22,545 14,46, 30,62,1,33 
Mie PUNTA 17-49-90 *47 95°57 ,5.37,21.53.15,45,29001, 3,55, 18 
Mie DATA S10 11.9 45.77 5957.59 525,55.15547, 31,65 

998 RESTORE 

B99 ° GO TO. 1 
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